N<-100	#入力次元
JRandSd<-1/N	#標準偏差
JRandM<-0	#平均
XRandSd<-1/sqrt(N)	#標準偏差
XRandM<-0	#平均
eta<- 1.0
M<-100

#set.seed(300)
b<-matrix(rnorm(N,JRandM,JRandSd),N,1)
#set.seed(500)
j1<-matrix(rnorm(N,JRandM,JRandSd),N,1)
#set.seed(800)
j2<-matrix(rnorm(N,JRandM,JRandSd),N,1)

R<-matrix(0,2*M+1,2)
L<-matrix(0,2*M+1,2)
tlearn <- function(){
	for (i in 1:M) {
		x<-rnorm(N,XRandM,XRandSd)
		ub<-x%*%b
		uj1<-x%*%j1
		j1<<-j1+eta*(ub-uj1)*x
		uj2<-x%*%j2
		j2<<-j2+eta*(ub-uj2)*x
		L[i,1]<<-sqrt(sum(j1*j1))
		L[i,2]<<-sqrt(sum(j2*j2))
		R[i,1]<<-sum(j1*b)/sqrt(sum(j1*j1)*sum(b*b))
		R[i,2]<<-sum(j2*b)/sqrt(sum(j2*j2)*sum(b*b))
		#R<-t(j1)%*%j2/(t(j1)%*%j1*t(j2)%*%j2)
		#R<-sum(j1*j2)/(sum(j1*j1)*sum(j2*j2))
	}
}
mlearn <- function(){
	for (i in (M+1):(2*M+1)) {
		x<-rnorm(N,XRandM,XRandSd)
		uj1<-x%*%j1
		uj2<-x%*%j2
		j1<<-j1+eta*(uj2-uj1)*x
		j2<<-j2+eta*(uj1-uj2)*x
		L[i,1]<<-sqrt(sum(j1*j1))
		L[i,2]<<-sqrt(sum(j2*j2))
		R[i,1]<<-sum(j1*j2)/sqrt(sum(j1*j1)*sum(j2*j2))
		R[i,2]<<-sum(j2*j1)/sqrt(sum(j2*j2)*sum(j1*j1))
		#R<-t(j1)%*%j2/(t(j1)%*%j1*t(j2)%*%j2)
		#R<-sum(j1*j2)/(sum(j1*j1)*sum(j2*j2))
	}
}
tlearn()
print(R,quote=F)
plot(L[,2],R[,2])
plot(1:(2*M+1),R[,1])
